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ABSTRACT 

The formation and evolution of a circumstellar disk in magnetized cloud cores 
is investigated from prestellar core stage until ~ 10"^ yr after protostar formation. 
In the circumstellar disk, fragmentation first occurs due to gravitational insta- 
bility in a magnetically inactive region, and substellar-mass objects appear. The 
substellar-mass objects lose their orbital angular momenta by gravitational in- 
teraction with the massive circumstellar disk and finally fall onto the protostar. 
After this fall, the circumstellar disk increases its mass by mass accretion and 
again induces fragmentation. The formation and falling of substellar-mass ob- 
jects are repeated in the circumstellar disk until the end of the main accretion 
phase. In this process, the mass of fragments remain small, because the circum- 
stellar disk loses its mass by fragmentation and subsequent falling of fragments 
before it becomes very massive. In addition, when fragments orbit near the 
protostar, they disturb the inner disk region and promote mass accretion onto 
the protostar. The orbital motion of substellar-mass objects clearly synchronizes 
with the time variation of the accretion luminosity of the protostar. Moreover, 
as the objects fall, the protostar shows a strong brightening for a short duration. 
The intermittent protostellar outflows are also driven by the circumstellar disk 
whose magnetic field lines are highly tangled owing to the orbital motion of frag- 
ments. The time-variable protostellar luminosity and intermittent outflows may 
be a clue for detecting planetary-mass objects in the circumstellar disk. 

Subject headings: accretion, accretion disks: ISM: clouds — stars: formation — 
stars: low-mass, brown dwarfs: planetary systems: planetary disks 
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1. Introduction 

The observation of many star-forming regions has shown that stars form in molecular 
cloud cores. A star is born after a long journey through the gravitational collapse of a cloud 
core. The gas begins to collapse in a dense part of the cloud (i.e., the cloud core) and increases 
its density and temperature as the cloud core collapses further. When the collapsing gas 
reaches a sufficiently high den sity (~ lO^^ c m"^), the collapse stops and a protostar having 



almost a Jovian mass is born (jLarsonlll969l ). Then, the protostar increases its mass through 



mass accretion and finally evolves into a main-sequence star (or main-sequence phase). 

The star formation process can be divided into two phases: the early protostellar collapse 
phase and later accretion phase of the star formation. The early phase, also called the gas 
collapsing phase, is defined as the period before protostar formation after gas collapse is 
initiated in the cloud core. Thus, during this phase, the gas continues to collapse. On the 
other hand, the later phase is defined as the period before the main-sequence phase after 
the gas collapse stops and a protostar appears in the collapsing cloud core. In particular, 
during the later phase, the period when the protostar significantly increases its mass by mass 
accretion is called "the main accretion phase." Since the star acquires almost all its mass 
in this main accretion phase, the final stellar mass is determined in this phase. Protostars 
during the main accretion phase are observed (or defined) as Class or I objects by spectral 
energy distribution (SED). The observation of star-forming regions also has shown that 
almost all Class or I objects drive a protostellar outfiow. Because a certain fraction 
of mass and angular momentum are ejected by the protostellar outfiows, they are closely 
related to the rate of mass accretion onto the protostar and final stellar mass. Moreover, 
the infalling gas with an angular momentum forms the circumstellar disk during the main 
accretion phase. The evolution of the circumstellar disk is related to the mass accretion and 
outfiow rates. In addition, the circumstellar disk is the site of planet formation. Therefore, 
understanding the evolution of the circumstellar disk during the main accretion phase is 
essential to understanding both the star and planet formation processes. 

Since the protostar and circumstellar disk during the main accretion phase are veiled by 
the infalling envelope, it is difficult to observe t hem. Recent l y, howe ver, the Spitzer Space 



Telescope (SST) has been unveiling this phase. lEnoch et al.l (l2009al ) discovered a massive 
circumstellar disk of ~ 1 Mq comparable to a central protostar around a Class object, 
indicating that (i) the disk already exists in the main a ccretion phase, and (ii) the disk mass 



is significantly larger than the theoretical prediction. lEvans et al.l ( l2009l ) observed several 
star forming regions and identified several Class objects using the SST. They showed that 
the bolometric luminosity of the objects is considerably dimmer than classical theoretical 
predictions, and has a dispersion over 2-3 orders of magnitude. They pointed out that recent 
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observations aggravate the "luminosity problem" (the p roblem that the accr etion luminosity 
of protostars is lower than theoretical predictions, see iKenyon et al.lll990[ ). and concluded 
that non-steady accretion is inevitably required to explain the observational results. The 
non-steady accretion may be attainable when the circumst ellar disk is sufficiently massive 
and causes the gravitational instability ( iDurisen et al.l 120071 ). 



Another r ecent development in observation of star and planet formatio n is direct images 



of exo-planets (IKalas et al. 



20081 : iMarois et al.ll2008l : iThalmann et al.ll2009r i , in which planets 



are locat ed at 10 AU from th e central star. However, in the framework of the core accretion 



scenario (iHayashi et al. Ill985l ). there is less likehood of planets forming in a region that is so 



remo te from the central star. Alternatively, the gravitational instability scenario (jCameron 
19781 ) may explain the formation of such planets. In summary, recent observations seem to 
indicate that the protostars have a considerably massive disk unlike what was previously 
believed. 

Although it seems that recent observational progress is unveiling problems for the early 
evolution stage of star formation, we cannot directly observe the circumstellar disk and pro- 
tostar in the main accretion phase because they are embedded in the dense infalling envelope. 
Thus, theoretical study is necessary to understand the properties of the circumstellar disk 
and protostar. In particular, multi-dimensional simulations are necessary to investigate the 
evolution of the circumstellar disk, protostellar outflow, and so on. 

The evolution of the collapsing gas cloud core from the protostellar core stage (i.e., 
the gas collapsing phase) unt il protostar formation has been well investigated using multi- 



dimensional simulations (e.g., Batelll998 



20071 : iBaneriee k Pudritzl l2006l : iMachida et al 



Tomisaka 



2002 



20061. 



Whitehouse k BatelboOG : IStamatellos et al 



2007 



2008a 



2009a). On the other 



hand, only a few studies have focused on the main accretion phase immediately following 
the prestellar core stage, because it is difficult to calculate long-term evolution in the main 
accretion phase with sufficient spatial resolution. In unmagnetized cloud cores, the forma- 
tion and evolution of the circumstell ar disk in the main ac cretion phase thro u gh th e gas 
collapsing phase were investigated by lWalch et al. fl2009al lb[). IVorobvov k Basul feoioh . and 
Machida et al.l (120101 ). They found that the circumstellar disk is considerably massive to in- 
duce fragmentation or the gravitational instability that is related to a non-steady accretion 
flow onto the protostar. 



In reality, however, since molecular clouds are strongly magnetized (jCrutcherlll999l ). the 
magnetic field may play ar i important role in the evolutio n of the circumstellar disk during 
the main accretion phase. IVorobyov k Basul (l2006l . l2007l ) investigated the evolution in two 
dimensions of the circumstellar disk in a magnetized cloud core and showed the non-steady 
accretion onto the central protostar. In three dimensions, the formation and evolution of the 
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circumstellar disk from prestellar core stage were investigated only by llnutsuka et al.l ( 120091 ). 
in which they showed fragmentation and possible planet formation in the magnetically in- 
active region of the circumstellar disk during the main accretion phase. They also indicated 
non-steady mass accretion onto the protostar owing to the gravitational instability of the 
circumstellar disk. However, this study only calculated the evolution of the circumstellar 
disk about ~ 1000 yr after protostar formation. 



In this study, in a setting similar to llnutsuka et al.l (|2009[ ). we investigate the evolution 
of the circumstellar disk for ~ 10'^ yr, which is ~ 1 times longer than the previous study. In 
addition to the model adopted by llnutsuka et al.l (120091 ). we newly calculate the evolution 
of the circumstellar disk formed in a relatively stable initial cloud core. In both models, we 
compare the mass accretion rate, properties of protostellar outflow, and the fragmentation 
condition. The structure of the paper is as follows. The framework of our models and the 
numerical method are given in §2. The numerical results are presented in §3. We discuss the 
fragmentation condition of the circumstellar disk and its implication for the planet formation 
in 54, and we summarize our results in 55. 



2. Model and Numerical Method 



2.1. Basic Equations 

To study the formation and evolution of a circumstellar disk in a magnetized molecular 
cloud core, we solve the three-dimensional resistive MHD equations, including self-gravity: 

(1) 



dv 1 
p— + p(v ■ V)v = - VP - —B X ( V X B) - pV0, 



dB 

~dt 



V X (v X B)+r]oY)V^B, 
V^(f) = AnGp, 



(2) 

(3) 
(4) 



where p, v, P, B, r^oo, and denote density, velocity, pressure, magnetic flux density, 
resistivity, and gravitational p otential, respectively. In addition, we adopted the hyperbolic 
divergence cleaning method of iDedner et al. I (l2002l ) to obtain divergence-tree magnetic field 
CV ■ B = 0). With this me t hod, no magnetic monopoles appear throughout calculation (see 



also Machida et al. 2005a: Matsumoto 



20071 ). 



To mimic the temperature evolution calculated bv iMasunaea fc fautsu 



adopt the piece- wise polytropic equation of state (see lVorobyov &: Basull2006 



ka (2000), we 



Machida et al. 
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2003) as 



P 



where Cs,o = 190 ms ^ and pc = 3.84 x 10 ^''g cm 



2/5 



[rir 



(5) 



10^" cm"'^). With equation (ED 



the gas behaves isothermally for n ^ 10 cm~ and adiabatically for n ^ 10 cm~ . For 
a reahstic evolution of the r nagnetic field in the c ircumstellar disk, we adopt the resistivity 
(^700) as the fiducial value in lMachida et al.l (120071 ). in which the Ohmic dissipation becomes 
effective for 10^^ cni ~^ ^ n <^ 10^^ cm~^ (for details, see Eq. [9] and [10], and Fig. 1 of 



Machida et al.ll2007f ). 



2.2. Initial Setting and Numerical Method 

As the initial state, we take a spherical cloud core with a critical Bonnor-Ebert (BE) 
density profile, in which a uniform density (pamb — 0.07pc) is adopted outside the sphere 
(r > Rc). For the BE density profile, we adopt the central density as ric = lO^cm"^ 
and isothermal temperature as T = 10 K. For these parameters, the critical BE radius is 
Rc = 4.7 X 10^ AU. The gravitational force is ignored outside the host cloud (r > Rc) to 
mimic a stationary interstellar medium. In addition, we prohibit the gas inflow at r = i?c to 
suppress mass input from the interstellar medium into the gravitationally collapsing cloud 
core. Thus, only the gas inside r < Rc collapses to form the circumstellar disk and protostar. 
Note that although the protostellar outflow driven by the circumstellar disk propagates into 
the region of r > Rc and disturbs the interstellar medium over time, we can safely calculate 
the mass accretion process onto the circumstellar disk because the outflow propagating into 
the interstellar medium does not affect the inner cloud core region (r < Rc). In this study, 
we call this initial spherical cloud with Bonnor-Ebert density profile 'the cloud core' that 
has a radius of Rc and mass of Md- 

Since the critical BE sphere is in an equilibrium state, we increase the density by a factor 
of / to promote the contraction, where / is the density enhancement factor that represents 
the stability of the initial cloud core. An initial cloud core with larger / is more unstable 
against gravity. In general, the stability of the cloud core is represented by a parameter 
ao (= E th/E^ra^r) that is the ratio of thermal (-Eth) to gravitational (-Egrav) energy. As 



shown in lMatsumoto fc Hanawal ( l2003l ). when the BE density profile is adopted, the density 



enhancement factor is related to the parameter ao as 

ao = — (6) 
We constructed two models with different / (= 1.68 and 2.8). In this paper, we call the 
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model having / = 1.68 [a = 0.5) "model AOS," and the model having / = 2.8 {a = 0.3) 
"model AOS." At the initial state, model AOS is more unstable than model AOS. 

In both models, the initial cloud core rotates rigidly around the z-axis in the region of 
r < Rc, while the uniform magnetic field parallel to the z-axis (or rotation axis) is adopted in 
the whole computational domain. In addition, we adopt that the angular velocity decreases 
in proportional to oc exp(— r^) outside the host cloud (r > Rc). Both models have the same 
angular velocity with Qq = 1.1 x 10~^^ s~^ and the same magnetic field, Bq = S7/iG. The mass 
inside r < R^ for each model is 0.8 (model AOS) and 1.3 Mq (model AOS). Model names 
and initial values are summarized in Table 1. The initially magnetized, rotating cloud core 
is also specified by the non-dimensional parameters (3q and 70, where f3o (= -Erot/|-E'grav|) and 
7o (= -Emag/I-Egravl) IS the ratio of the rotational and magnetic energy to the gravitational 
energy inside the initial cloud core. Model AOS has /3o = S x 10^'^ and 70 = O.OS, while 
model AOS /3o = 6 x 10^^ and 70 = 0.07. In addition, the mass-to-flux ratio M/$ is also 
used to specify the initial cloud, where M and $ are the mass and magnetic flux of the initial 
cloud core. There exists a critical value of M/$ below which a clo ud is supported against 



the gr avity by the magnetic field. For a cloud with uniform density, iMouschovias fc Spitzer 



(119761 ) derived a critical mass-to-fiux ratio 



where the constant ( = O.SS. We define the mass-to-fiux ratio normalized by the critical 
value as A. In our setting, model AOS has A = 9, while model AOS has A = S.6. 

We add m = 2-mode non-axisymmetric density perturbation to the initial cloud core. 
Then, the density profile of the cloud core is described as 



p(r) 



PBE(r) (1 + 6p) f for r < Rc, , . 

PBE{Rc){l + Sp)f for r>R,, ^ ^ 

where Pbe('") is the density distribution of the critical BE sphere, and Sp is the axisymmetric 
density perturbation. For the m = 2-mode, we chose 

5p = A^{r/R,fcos2(l), (9) 

where (=0.01) represents the amphtude of the perturbation. The radial dependence is 
chosen so that the density perturbation remains regular at the origin (r = 0) at one time-step 
after the initial stage. This perturbation ensures that the center of gravity is always located 
at the origin. 

In the collapsing cloud core, we assume the protostar formation to occur when the 
number density exceeds n > 10^^ cm~^ at the center. To model the protostar, we adopt a sink 
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around the center of the computational domain. In the region r < rgink = 1 AU, gas having 



a number density of n > 10^^ cm is removed from the c omputational domain and ad ded to 



-3 ; 



the protostar as gravity in each timestep (for details, see lMachida et al.ll2009al . l2010| ). This 
treatment of the sink makes it possible to calculate the evolution of the collapsing cloud 
core and circumstellar disk for a longer duration. In addition, inside the sink, the magnetic 
flux is removed by Ohmic dissipation, because s uch a region has the magnetic Reynolds Re 
number exceeding unity Re > 1 (for details, see lMachida et al.l 120071 1. 



see 



To calculate on a large spatial scale, the nested grid method is adopted (for details. 



Machida et al. ll2005al Jbl). Each level of a rectangular grid has the same number of cells 
(128 X 128 X 32). The calculation is first performed with five grid levels (/ = 1 — 5). The 
box size of the coarsest grid Z = 1 is chosen to be 2^i?c- Thus, grid of Z = 1 has a box 
size of ~ 1.5 X 10^ AU . A new finer grid is generated before the Jeans condition is violated 
( iTruelove et al. 1119971 ). The maximum level of grids is /max = 12 that has a box size of 74 AU 
and a cell width of 0.58 AU. 

We adopted a fixed boundary condition on the outermost grid boundary. The uniform 
density (pamb) and magnetic field {B^ = 0, By = 0, B^ = Bq) and zero fluid velocity 
(vx = Vy = Vz = 0) are imposed on / = 1 grid boundary at each timestep. Such boundary 
condition hardly affects the evolution of the cloud core, because the gravitationally collapsing 
cloud core with a radius of 4.7 x 10^ AU (r = Rc) is embedded in a large simulation box with 
a size of ~ 3 x 10^ AU, inside which the static interstellar medium is fulfilled in the region 
of r > Rc- In our setting, the Alfven speed in the interstellar medium (r > tbe) for model 
A03 is va = 0.15kms~^. Thus, by reaching the Alfven wave to the computational boundary 
from the center of the cloud core, it takes 4.7 x 10® yr. Since we stopped our calculation in 
~ 10^ yr after the calculation begins, the Alfven wave generated at the center of the cloud 
core (or the computational boundary) never reaches the computational boundary (or the 
center of the cloud core). 



3. Results 



We calculated the formation and evolution of the circumstellar disk in magnetized col- 
lapsing cloud cores. First, we simply outline the formation and evolution of the circumstellar 
disk. Before protostar f ormation, the first (adiabatic) core forms in the collapsing cloud cor e 
with a scale of ~ 10 AU (lLarsonlll969l : iMasunaga fc Inutsuka Il2000l : ISaigo fc Tomisakall2006l ) . 
At its formation, the first core has a thick disk-like structure and is mainly supported by the 
thermal pressure gradient force. After protostar formation, the first co re becomes thin be- 
cause the centrifugal force dominates the thermal pressure gradient force (IWalch et al.ll2009al : 
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Machida et al.ll2010l ). During the main accretion p hase, the f irst core grows and ex tends to 



a large extent (~ 100 - 1000 AU). Ilnutsuka et al.l fcOOQf ) and Machida et all fcOlOf ) pointed 
out that the circun astellar disk origi nates from the first core, and is formed before protostar 
formation (see also Bate IQQsl . 2010). Their calculations showed that the first core gradually 
grows to become the circumstellar disk in the main accretion phase. They also pointed out 
that, reflecting the thermal history of the collapsing cloud core, the circumstellar disk is in- 
evitably more massive than the protostar in the early part of the main accretion phase. The 
formation of the first core (or the circumstellar disk) before protostar formation has been 



inves tigated in many past studies (see review by iBodenheimer et al. Il2000l : iGoodwin et al. 



20071 ). Thus, although we began the calculation from the prestellar core stage, we mainly 



focus on the evolution of the circumstellar disk in the main accretion phase in this paper. 



3.1. Recurrent Fragmentation and Planet Formation 



As described in ^ we constructed two models (models AOS and AOS) with different 
initial stabilities (i.e., different ao). In the two models, the protostar forms t ~ 2.50 x 10^ yr 
(model AOS) and t ^ 5.48 x 10^ yr (model A05), respectively, after the calculation begins. 
Thus, in the collapsing cloud core, the protostar for model AOS forms for a shorter duration 
than that for model A05. This is because the initial cloud core for model AOS is more unstable 



than that for model A05 and begins to collapse acco rding to the self-simila r solution (ILarson 
19691 ) in a shorter time after the calculation begins (IMachida et al.ll2008al ). In both models, 
we calculated the evolution of the circumstellar disk for ~ 10^ yr after protostar formation. 

First, we describe the evolution of the circumstellar disk for model AOS that shows 
the recurrent fragmentation and formation of substellar-mass object in the main accretion 
phase. Figure [1] shows the evolution of the circumstellar disk after protostar formation for 
model AOS, in which the density distributions around the protostar on the equatorial plane 
are plotted. The box sizes in the upper panels are ~ 80 AU, while those in the middle 
and lower panels are ~ S20AU. Figure [T]a shows the structure of the circumstellar disk 
tc — 450 yr after protostar formation, in which the red region {n ^ 10^^ cm~^) corresponds 
to the circumstellar disk. Note that we describe the elapsed time after protostar formation 
with a notation of 'tc,' which differs from the elapsed time 't' after the calculation begins 
(or the cloud begins to collapse). The size of the circumstellar disk increases with time and 
extends up to ~285 AU by the end of the calculation. In the circumstellar disk, fragmentation 
o ccurs due to gravitational instability tg — 6S0 yr after protostar formation (for details, see 



Inutsuka et al.ll2009l : iMachida et al.ll2010l ). Two ambiguous clumps appear in Figure [Tf> and 
c, while two clear fragments are seen in Figure [T](i. Fragmentation occurs ~ 6.8 AU away 
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from the protostar. Then, fragments orbit around the protostar with an orbital separation 
of ~ 2 — 38 AU (Figs. [T]c-e) and finally fall onto the protostar at — 3760 yr (Figs. [T]e and f). 
In this section, we call the fragment appeared in the circumstellar disk the planet. Note that 
in reality, since some fragments that appeared in this calculation exceeded the deuterium- 
burning limit (about 13 Jupiter mass), they are not planets in the complete sense. Note 
also that we define, for convenience, the objects born in the circumstellar disk as planets to 
characterize them as a whole: we discuss their mass in §4.21 

Figure [2] shows the structure of the circumstellar disk and protostellar outflow at almost 



the collapsing cloud core has been studied by many authors (e.g., 


Matsumoto & Tomisaka 


2004; 


Machida et al. 


2004; 


Machida et al. 


2005a; 


Hennebelle & Fromang 


2008a 


). Even in 



the main accretion phase (i.e., even after protostar for mation), the circumstellar disk that 



originates from the first co r e conti nues to drive outflow ( lTomisakall2002l ; iBanerjee fc Pudritz 



20061 : iMachida et al.l 120071 . |2009aJ). In this model, the outflow driven by the circumstellar 
disk extends up to ~516AU by the end of the calculation. As seen in Figure |2l the outflow 
is driven only by the outer disk region, while the outflow is not driven by the inner region of 
the circumstellar disk where fragmentation occurs and planets appear. This is because the 
inner disk region has a density of n >, 10^^ — 10^^ cm~^ a nd th e magnetic field is dis sipated 



by Ohmic dissipation. As shown in iNakano et al.l ( 120021 ) and iMachida et al.l ( 120071 ). when 
the number density exceeds n ^ 10^^ — 10^^ cm~^, Ohmic dissipation becomes effective owing 
to the extremely low ionization degree. Note that the ambipolar diffusion can be effective in 
the range of n < 10^^ — 10^^ cm~^ (see, §4.3p . Fragmentation is suppressed by the magnetic 
field because magnetic effects such as magnetic bra king and outflows effectively transfer 



the angular momentum that promote s fragmentation (IMachida et al. Il2005bl : IMachida et al. 



2008al : iHennebelle fc Teyssied l2008bl ). Thus, a weaker field promotes fragmentation in the 
circumstellar disk. As a result, the inner disk region has a considerable weak magnetic field 
and fragmentation occurs without outflow, while the outer disk region has a strong field and 
no fragmentation occurs with outflow. 

In Figure [2], yellow lines are magnetic field lines. In the figure, only magnetic field lines 
that threaded planets are plotted to stress the planet's spin or orbital motion. The magnetic 
field lines just above each planet are strongly twisted by the spin motion of the planet. In 
the midair (or on a large scale), the magnetic field lines originating from each planet are 
twisted together, reflecting the orbital motion of the planets. Although we only plotted 
magnetic field lines near the planet to understand the planet's motion, the magnetic field 
(or Lorentz force) in this region (or inner disk region) is extremely weak (/3p ^ 10, where 
Pp is the plasma beta and different from initial ratio of the rotational to the gravitational 
energy, Pq) because Ohmic dissipation is effective. Thus, near the planets, the magnetic field 
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lines are passively moved, and they hardly affect the dynamical evolution of both the planet 
and circumstellar disk. However, planets perturb the outer disk region that is connected 
to magnetic field lines driving outflow (i.e., the magnetic field lines in the outer disk region 
where Ohmic dissipation is not effective). As a result, a highly time- variable outflow appears 
in this model. We will discuss an intermittent outflow in §3.71 

Figure [T]e shows the structure of the circumstellar disk just before the planets fall onto 
the protostar, while Figure [TJf shows the disk just after the fall. After the first-generation 
planets disappear, fragmentation occurs again in the circumstellar disk and the second- 
generation planets appear tc — 3990 yr after protostar formation as shown in Figure [T]g. 
However, these planets also fall onto the protostar ~ 680 yr after their formation. The 
formation and falling of the planets are repeated several times in the circumstellar disk. 
By the end of the calculation, 18 planets appeared and 16 planets fell onto the protostar. 
Figure [3] shows eighth- and ninth-generation planets. The eighth-generation planets fell onto 
the protostar, while the ninth-generation planets survived until the end of the calculation. 
The orbital radius r^, time tj at planet formation, planet's falling epoch te, lifetime tg — ti, 
and mass at its formation Mp\ are summarized in Table 2. 



3.2. Orbital Trajectories of Planets 

Figure H] shows the orbital evolution of planets against time after protostar formation 
tc- When multiple planets appear as seen in Figures [1] and [3l only the orbital radius of the 
most massive planet is plotted. The figure indicates that the first-generation planets form at 
~ 6.8 AU from the protostar at — 795 yr, and shrink their orbit and fall onto the protostar 
after ~ 10 orbital rotations. The second-generation planets appear at ~ 46 AU and fall 
onto the protostar after ~ 2 — 3 orbital rotations. This formation and falling of planets is 
repeated, until the ninth-generation planet is formed. 

Figure shows the trajectory for planets for each generation. The figure indicates that 
planets shrink their orbits with a certain amount of eccentricity. Figures S] and show 
that planets are born at ~ 6 — 50 AU from the protostar and orbit for ~ 500 — 1000 yr (or 
2 — 10 orbital rotations) around the protostar. In addition, these figures show that the first- 
generation planets are formed at a relatively small orbital radius (~ 6.5 AU) and orbit for 
~ 10 times around the protostar, while other generations of planets are formed at a larger 
orbital radii (~ 20 — 50 AU) and orbits only for 2 — 3 times. This is because the circumstellar 
disk increases its size and mass with time. As shown in Figure [H the circumstellar disk has 
a size of ~ 20 — 30 AU when the first-generation planet appears, while the disk size reaches 
~ 200 AU just before the second-generation planet appears. As a result, the magnetically 



inactive zone also expands as the disk mass and radius increase, and fragmentation and planet 
formation can occur in such an outer disk region. The mass evolution of the circumstellar 
disk is shown in §3.4^ and the fragmentation condition is discussed in §3.61 



3.3. Non-Fragmentation Model 

Unlike model A03, model AOS shows no fragmentation in the circumstellar disk by the 
end of the calculation. Figure |6] shows the density distribution on the y = (left) and z = 
(right) planes = 3781 yr after protostar formation for model AOS. As in model A03, model 
AOS also shows the outflow driven by the circumstellar disk. The outflow appears ~ 700 yr 
before protostar formation and reaches ~ 4400 AU at the end of the calculation. In the right 
panel, the ring-like density gap is seen in the range of 40 AU <^ r ^ 80 AU. In this range, 
the gas is largely blown away from the circumstellar disk by the outflow. 

In the right panel, a non-axisymmetric structure is seen in the proximity of the protostar 
(r ^ 40 AU), while nearly axisymmetric structures are seen in the outer disk region (r ^ 
40 AU). The magnetic field is well-coupled with neutrals in the outer disk region because the 
gas density is lower than n <^ 10^^ cm~^ where Ohmic dissipation is not effective. Thus, the 
angular momentum is effectively transferred by magnetic braking in the outer disk region, 
and the gas in such a region can fall into the inner disk region. On the other hand, the 
magnetic field is significantly dissipated by the Ohmic dissipation in the inner disk region 
(r <^ 40 AU), because the density in such a region exceeds n ^ 10^^ cm~^. Thus, in the inner 
disk region, the angular momentum transfer by the magnetic field is not effective, and the 
excess angular momentum and ac cumulated mass cau s e the gravitational ins tability to induce 



the non-axisymmetric structure ( IWalch et al.ll2009ai iMachida et al.l 120101 ). However, since 



the increase rate of the disk mass for model A03 is considerably smaller than that for model 
AOS (see, §3.4p . the circumstellar disk for model AOS is not sufficiently massive to induce 
fragmentation (see, §3.6p . For model AOS, fragmentation and formation of substellar-mass 
object do not occur until the end of the calculation. As a result, since the circumstellar disk 
is not disturbed by the fragments, the protostellar and circumstellar-disk mass moderately 
increase, unlike in the case of model A03 (see, §3.4p . 



3.4. Masses of Protostar and Circumstellar Disk 



Figure [7] shows the evolution of protostellar mass against time after protostar forma- 
tion. The protostellar mass Mpg is derived as the mass falling into the sink. At its formation. 



- 12 - 



the pr otostar has a mass of Mps ~ lO'^Mp, as expected in spherically symmetric calcula- 
tions (lLarsonlll969l : iMasunaga fc Inutsuka Il2000l ). In the main accretion phase, protostars 
increases their mass by gas accretion and reach Mpg = 0.59 (model AOS) and 0.15 
(model A05), respectively, at the end of the calculation. In Figure [3, the step-like increase of 
the protostellar mass for model AOS is caused by the falling of the planet. Since the planet 
and its envelope have a total mass of ~ 6 — 27Mjup (see. Table 2), the protostar drastically 
increases its mass at the moment of the falling. In both models, the protostars have almost 
the same mass for <^ SOOO yr, while the protostar for model AOS is more massive than that 
for model A05 for tc ^ SOOO yr. The difference in protostellar mass is caused by the difference 
in the rate of mass accretion onto the protostar. The different mass accretion rates are ow- 
ing to different efficiencies of the angular momentum transfer in the circumstellar disk. The 
higher rate of mass accretion onto the protostar is realized in a more massive circumstellar 
disk, because a massive disk tends to show gravitational instability that effectively transfers 
the angular momentum outward, and increases the rate of mass accretion onto the protostar. 

Figure [8] shows the time evolution of th e circumstellar-disk mass that is estimated in 
the same manner as in iMachida et al.l ( l2010l ) . The figure indicates that the disk mass for 
model AOS is significantly greater than that for model A05 at the epoch just after protostar 
formation. The circumstellar disk for model A05 gradually increases its mass until the end 
of the calculation and reaches Mdisk ~ 0.042 for ~ 10^ yr. The average increase rate of 
disk mass for model A05 is Mdisk ~ 4.2 x 10~^Mq yr^^. On the other hand, the circumstellar 
disk mass for model AOS reaches M^isk ~ 0.15 Mq for ~ 10^ yr with an average increase rate 
of Mdisk ~ 1-5 X lO~^M0yr~^. Thus, considering that the initial stability of the cloud core 
is closely related to the mass accretion rate (see, §4.21) . the increase rate of the disk mass 
for model AOS is about one order of magnitude higher than that for model A05. Note that 
it is difficult to accurately estimate the total mass of accreting gas onto the circumstellar 
disk because a fraction of the circumstellar-disk mass falls onto the protostar with time. For 
model AOS, the disk mass shows a high time variability that is caused by the time variable 
mass accretion onto the protostar (see, §S.5p : the spike-like decrease of the disk mass for 
model AOS corresponds to the epoch at which the planet falls onto the protostar. 



3.5. Rate of Mass Accretion onto Protostar 

Since the gas accretes onto the protostar mainly through the circumstellar disk in the 
main accretion phase, the mass accretion rate depends on disk properties such as disk mass, 
spiral structure, and existence of fragments. Figure [9] shows the mass accretion rates for 
~ 10^ yr after protostar formation. We estimated the mass falling into the sink every one 
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year and defined it as the rate of mass accretion onto tlie protostar, Mpg. Tlius, tfie time 
resolution of tfie accretion rate in Figure |9] is one year, while the time-step of calculations is 
At ^ 0.001 — 0.01 yr. Therefore, we integrated the falling mass every ~ 100 — 1000 time-steps. 

In Figure [H the mass accretion rate for model AOS gradually decreases, oscillating 
between M ~ 10~^ and M ~ 1O~^M0 yr~^, while model AOS shows intermittent gas accretion 
with multiple peaks. For model A05, it is considered that the spiral structure in the inner 
disk region induces the non-steady accretion flow. The spiral structure can transfer the 
angular momentum outward and move the gas inward. The range of the accretion rate in 
this model almost corresponds to the theoretical prediction that is expressed as 

M = ~ 1.9 X lO-VmM0yr-\ (10) 



where /m is a numerical factor (e.g., /m = 0.975 for IShulll977l . /m = 46.9 for iHunted 119771 ) . 
and Cs = 200 ms~^ is adopted. The gas falls onto the circumstellar disk first, and the rate 
of gas accretion onto the circumstellar disk is ~ 4 x 1O~^M0 yr~^ on average, as shown in 
Figure [H In summary, the rate of mass accretion onto the protostar almost corresponds to 
the mass accretion onto the circumstellar disk for model AOS. 

On the other hand, model A03 shows intermittent mass accretion caused by the large- 
scale perturbation of the disk and falling of planets. As shown in Figure HI the first-generation 
planets form at r ~ 6.8 AU and orbit in the range of 2 to 40 AU until they fall onto the planet. 
They disturb the circumstellar disk and transfer the angular momentum outward. Thus, the 
mass accretion rate is linked to the orbital motion of the planet especially for tc <^ 3000 yr, i.e., 
before the first-generation planets fall onto the protostar. After the first-generation planets 
disappear, planets of subsequent generations (2nd to 9th) appear at r 20 — 50 AU from 
the protostar. They orbit in the range of from 5 to 50 AU and approach the protostar within 
r < 5 AU only before falling onto the protostar. Thus, just before falling, these planets can 
disturb the inner disk region and promote mass accretion onto the protostar. The peak mass 
accretion rate for tc ^ 3000 yr corresponds to the epoch in which a planet approaches and 
falls onto the protostar. In addition, since not only the gas in the circumstellar disk but also 
the planet itself falls onto the protostar, the mass accretion rate at its peak for a model with 
planets is much larger than the theoretical prediction. In Figure [9|, there are eight peaks of 
the mass accretion rate for tc ^ 3000 yr, corresponding to epochs of falling of planets. In such 
an epoch, the protostar has a very high mass accretion rate of Mps ^ 10^^ — lO^^MQyr^^, 
because a planet with mass of ~ 10^^ — 10^^ Mq falls onto the protostar in a short duration. 
Note that this high accretion rate may be smoothed over time if calculated with higher 
spatial resolution, while it is considered that the total mass falling onto the protostar does 
not change. 
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3.6. Disk Properties for Models with and without Fragmentation 

Figure [TO] show s the radial d istribution of the surface density o (Fig. [TOja). Toomre Q 
parameter (Fig. [TOb. lToomrelll964j ) . plasma /3p (Fig.lTOb). and ratio of the azimuthal to Kepler 
velocity f</,/fkep (Fig. ITOb). These quantities are azimuthally averaged. The azimuthal- 
averaged surface density o [= cr(r)] is derived as 

air) = drt-sirO)r^de ^ ^^^^ 
do 

where, 

(Tsir,9) = j pir,e,z)dz, (12) 

and 

dS = dr / r^de, (13) 
Jo 

where Tc is the cylindrical radius, and dr = ds{l) (the cell width of Zth-level grid) is adopted. 
To determine Zcri in equation ( fT2|) . we estimated the disk height by eyes and 2;cri=10AU is 
adopted. We confirmed that the surface density a hardly depend on z^ri in the range of 
Zcri > 5 AU. In Figures [TOfc and c, Q parameter and plasma f3p are integrated along 2;— axis 
and azimuthally averaged. To estimate Q parameter, firstly, we integrated the mass-weighted 
sound speed Cs and epicyclic frequency k along z— axis as 

u.»)^^^^^=^^^^44?^^^. (14) 

and 

,^(,,),r'■"''M■;)M''■''--)''^ (15) 



Then, using equations f[T2|) . f[T4j) and f[T5|) . we derived the azimuthal-averaged Q parameter 



where 



as 

Q(r) ^ ^'-tQ^^X^.O ^ ^^^^ 

Q(r.«) = H4i^lfM. (17) 

The azimuthal-averaged plasma /3p(= SticIp/B"^) is also derived in a similar way to Q pa- 
rameter as 

fUr) = "■tPM-.l^ ^ (13) 
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where 

/?P>.9) = £^^^^^^^44f^^^. (19) 

In addition, the ratio of azimuthal-averaged f ^ to the Kepler velocity on the equatorial plane 
is plotted in Figure [TOk i. The black lines in Figure [10] are the radial quantities at the epoch 
just before the first fragmentation {tc = 445 yr) for model AOS, while blue lines are those 
at the epoch just before the second fragmentation after the first fragments falling onto the 
protostar (tc = 3947 yr) for model AOS. The red-dotted lines are quantities at almost the 
same epoch {tc = 3781 yr) as blue lines but for model A05. 

The black lines in Figure [TUla and b show that, for r < 10 AU, the surface density has 
a maximum and Q parameter becomes Q < 1. Just after this epoch, fragmentation occurs 
at r = 6.8 AU where Q < 1. The black lines in Figures [TOb and d show that plasma /3p 
suddenly rises at r ~ 10 AU inside which the flow moves azimuthally with near the Keplerian 
velocity. Since the gas (surface) density is sufficiently high in the region of r < 10 AU, the 
Ohmic dissipation becomes effective and magnetic field considerably dissipates. In general, 
fragmentation hardly occurs when the neutral gas is strongly coup led with the magnetic 
field (or ions), because the m agnetic field suppresses fragmentation (IHennebelle &: Teyssier 



2008bl : iMachida et al.ll2008bl ). Instead, when Q < 1, fragmentation is easy to occur in the 



magnetically inactive region where the Ohmic dissipation weakens the magnetic field. In 
summary, at this epoch, fragmentation can occur in the region of r < 10 AU, because Q < 1 
is realized owing to the higher surface density and very weak magnetic field (or higher /3p) 
does not suppress fragmentation. 

In the region of r < 10 AU, the surface density at = 3947 yr (blue line) is lower 
than that at tc = 445 yr (black line). This is because, in the inner disk region, the angular 
momentum is transferred by the orbital motion of the first-generation planets and the gas 
falls onto the protostar. In addition, planets c ontrib ute to the gap formation around the 
central protostar as seen in IVorobyov fc Basij (l2010l ). Thus, at tc = 3947 yr, the surface 
density near the protostar lowers and Q parameter exceeds Q > 1 for r ^ 30 AU (blue line 
of Fig. [TUb). Since the lower (surface) density makes the Ohmic dissipation ineffective, lower 
plasma /3p (/3p ~ 2 — 4) is realized in this region (r ^ 40 AU). Therefore, the inner disk region 
becomes magnetically active and the angular momentum is transferred also by the magnetic 
braking. Figure [TOH shows that azimuthal velocity to the Kepler velocity at tc = 3947 yr 
is lower than that at tc = 445 yr. The lower surface density and higher value of plasma 
(3p are not adequate to induce fragmentation. On the other hand, the surface density and 
plasma /3p at tc = 3947 yr have a local peak at r ~ 30 — 40 AU (blue lines of Figs. [TOl a and 
c). This is because, in this region, the radially inward motion of the flow is suppressed by 
the first-generation planet and the gas is accumulated. In such high-density gas region, the 
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Ohmic dissipation becomes effective owing to the higher surface density. Thus, the magnetic 
braking becomes ineffective and azimuthal velocity approaches the Kepler speed (blue line of 
Fig. [TUlrf) ■ Therefore, the centrifugal barrier further amplifies the surface density. Finally, Q 
parameter becomes Q < 1 in the range of 20 AU ^ r ^ 60 AU (Fig. [TUb). and fragmentation 
occurs to form the second-generation planets at r ~ 46 AU as seen in Figure [T]g. 

In Figures [TOla and 6, model AOS (red-dotted line) has a relatively lower surface density 
because this model has a lower mass accretion rate. Therefore, Q parameter never reach 
Q ^ 1, and fragmentation does not occur until the end of the calculation for model A05. 
In addition, the lower surface density increases the magnetically active area (Fig. [TOb ) and 
stronger outflow appears as seen in Figure [61 Thus, the excess angular momentum is mainly 
transferred by the protostellar outflow and magnetic braking for model AOS, while the excess 
angular momentum contributes to fragmentation for model A03. 

In summary, the initial stability (or mass accretion rate) is related not only the surface 
density (or mass of the circumstellar disk) but also the magnetic activity in the circumstellar 
disk. The higher accretion rate realizes the higher surface density (i.e., lower Q) and has 
a large magnetically inactive area where the angular momentum is not transferred by the 
magnetic effects. Both effects promote fragmentation and planet formation in the circum- 
stellar disk. Instead, since the lower accretion rate has a lower surface density and large 
magnetically active area, fragmentation and planet formation tend to be suppressed. 

3.7. Protostellar Outflow Driven by Circumstellar Disk 

The condition of the circumstellar disk also affects protostellar outflow. Figure [TT] shows 
the mass outflowing from the protostellar system (i.e., protostar + circumstellar disk system) 
for both models. To measure the outflowing mass, we integrated the gas having Vr > c^.o 
inside the initial radius of cloud core (r < i?c), where Vr is the radial velocity. The outflowing 
mass for model AOS gradually increases for ^ 6000 yr and decreases for tc ^ 6000 yr. 
Since the freefall timescale of the initial cloud core is ~ 1.8 x 10^ yr, it is expected that the 
outflow weakens after passing the peak and disappears on a freefall timescale. In addition, 
the outflow extends far beyond the initial radius of cloud core for tc ^ SOOOyr. Thus, we 
underestimated the outflowing mass for tc ;^ SOOOyr because we integrated the outflowing 
mass only inside r < Rc. To conclusively estimate the outflowing mass and total ejected 
mass from the initial cloud core, more long-term calculations are necessary. However, since 
such calculations are beyond the scope of this study, we do not comment on them any more. 
At its peak, the outflowing mass is about O.OSMq which is comparable to the protostellar 
mass Mps ~ 0.1 (see Fig. [7]) at the same epoch. Thus, a large fraction of the infalling 
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gas is blown away from the circumstellar disk by the protostellar outflow for model AOS. 

Unlike model AOS, in Figure [TTl the outflow mass for AOS shows high time variability, 
indicating that outflow is unsteadily driven by the circumstellar disk. In summary. Figure ITT] 
indicates that steady outflow is driven by the less massive circumstellar disk without frag- 
ments (model AOS), whereas intermittent outflow is driven by the more massive circumstellar 
disk with fragments (model AOS). Note that Figure [TT] shows only the outflowing mass hav- 
ing the high- velocity component of f > o, and ignores the low- velocity component of the 
outflow. 

Figure [12] shows the ratio of the outflowing to infalling mass. To investigate the outflow 
efficiency, we calculated the outflowing/infalling mass rate on a / = 9 grid surface, which 
covers the whole circumste llar disk w i th a g rid size of ~ 600 AU. The box size of / = 9 grid 



is ~ 700 AU. According to iTomisakal ( 120021 ) . we estimated the outflow Mout and inflow Mjj 

masses as 



Mout = / pmax[v ■ n, 0] dS, (20) 

J boundary of 1=9 



Min = / pmaxf-y ■ -n, 0] ^5*, (21) 

J boundary of 1=9 



and 



where n is the unit vector outwardly normal to the surface of the / = 9 grid. In Figure | 
model AOS keeps an almost constant ratio, Mout/^in ~ O.S, until the end of the calculation. 
Thus, about S0% of the infalling mass is ejected by the outflow. In this model, the outflow 
is driven by the outer disk region where the axisymmetric structure is kept, while no outflow 
appears in the inner disk region where a considerably weak magnetic field is realized by 
Ohmic dissipation and a non-axisymmetric structure appears, as seen in Figure |6l In the 
outer disk region, the disk mass is blown away by the outflow, and effective magnetic braking 
moves the gas inward. Therefore, such a region does not maintain mass sufficient to induce 
the gravitational instability. As a result, the steady disk region drives the steady flow for 
model AOS. 

On the other hand, the outer disk region for model AOS does not remain stable, because 
the increase rate of the disk mass (or the mass accretion rate onto the disk) for model AOS 
is much larger than that for model AOS. Although the magnetic braking and mass ejection 
by the outflow are effective in the outer disk region for model AOS also, they cannot transfer 
sufficient mass to suppress the gravitational instability and subsequent fragmentation owing 
to the high mass accretion onto the disk. Thus, multiple planets appear in the disk. Oscil- 
lation of the infall-outflow mass ratio for model AOS is shown in Figure [T21 Comparison of 
Figure m with Figure [T2] indicates that this oscillation is related to the orbital motion of the 
planet. It is considered that the magnetic field lines and surface density of the circumstellar 



- 18 - 



disk are disturbed by the planet's orbital motion, which causes highly time- variable outflow, 
as shown in Figure [T2j The time- variable protostellar outflow may be a clue for detecting a 
planet in a very early phase of star formation. 



4. Discussion 

4.1. Can We Detect Planets or Substellar-Mass Objects in the Main 

Accretion Phase? 

So far, it is considered that planets form in relatively quiet circumstellar disks after 
the gas accretion onto the circumstellar disk has almost ended. In this study, however, we 
showed the possibility of the formation of planets or substellar-mass objects by gravitational 
instability in a vigorous circumstellar disk during the main accretion phase, in which the 
circumstellar disk continues to increase its mass by gas accretion from the infalling envelope. 
When planets or substellar-mass objects appear in the circumstellar disk, they can affect 
the rate of accretion onto the protostar. As described in §3.5[ we may verify their existence 
by observing the time-variability of the rate of accretion onto the central protostar. As 
seen in Figure [HI the protostar shows a relatively low time- variability of the accretion rate 
when neither planets nor substellar-mass objects exist in the circumstellar disk, while it 
shows a high time- variability of the accretion rate when fragmentation occurs and planets or 
substellar-mass objects appear. Figure [13] shows the time variation of accretion luminosity 
of the protostar, which is deflned as 

L^^gM^, (22) 

where we adopted the protostellar mass Mpg and mass accretion rate Mpg as values derived 
from calculations at each time and used a constant protostellar radius rpg = 2 Rq for simplic- 
ity. The flgure indicates that the protostar for model AOS has an accretion luminosity in the 
range of 0.01 ^ Lps/ Lq ^ 1 for ~ 10^ yr after protostar formation. The accretion luminosity 
for this model gradually increases for a longer timescale of ~ 10^ yr, while it oscillates with a 
period of ~ 100 yr, which corresponds to the Kepler timescale in the region of a developing 
spiral structure. The Kepler timescale Tk can be expressed as 

As shown in Figure O the non-axisymmetric structure appears a t r ^ 10 — 40 AU, which 



corresponds to the magnetic non- active zone ( iMachida et al.l 120071 ). Thus, it is natural that 
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the protostar shows time-variable accretion luminosity with a period of ~ 100 yr, because 
the angular momentum transfer that is closely related to the rate of mass accretion onto 
the protostar is caused by the spiral structure in the circumstellar disk. The variation in 
amplitude of the accretion luminosity for a short duration (~ 100 yr) is about one order of 
magnitude. 

On the other hand, the accretion luminosity for model AOS is qualitatively and quanti- 
tatively different from that for model AOS. The protostar for model AOS repeatedly shows 
a strong brightening. In this model, the accretion luminosity suddenly increases to reach 
~ 10 — 10^ Lq when planets or substellar-mass objects approach the protostar or fall onto 
the protostar, while the luminosity is lower, ^ 0.1 — 0.01 Lq when planets or substellar-mass 
objects orbit far from the protostar. Thus, the strong brightness may be a clue for detecting 
planets or sub-stellar-mass objects. 

In addition, when the planets or substellar-mass objects orbit in the inner disk region of 
^ 10 — 20 AU, their orbital motion is linked to the accretion rate (or accretion luminosity). 
The accretion rate and orbital radius for t^. = 1000 — 2800 yr for model AOS are plotted 
in Figure [HJ In this time period, the first-generation planet orbits around the protostar 
for ~ 10 times. The figure indicates that the orbital motion of planets or substellar-mass 
objects clearly synchronizes with the mass accretion rate. The orbiting objects disturb the 
circumstellar disk and promote mass accretion onto the protostar. Since the circumstellar 
disk has a higher density and shorter timescale in the proximity of the protostar, a high mass 
accretion rate is realized when the planets or substellar-mass objects approach the protostar. 
As a result, the mass accretion rate is closely related to the orbital motion of planets and 
substellar-mass objects. We may estimate and verify the orbital motion of planets in the 
circumstellar disk by observing the time-variable mass accretion. Note that we ignored the 
effect of the infalling envelope when the accretion luminosity is estimated in Figure [T31 
Note also that, to quantitatively estimate the accretion luminosity, we have to consider the 
infalling envelope that weakens the luminosity from the protostar. 



4.2. Mass Accretion Rate onto Disk and Masses of Fragments 

In this study, we showed that fragmentation and subsequent formation of planet and 
sub-stellar mass object can occur in a massive circumstellar disk, while no fragmentation 
occurs in a relatively less massive disk. The mass of the circumstellar disk or the mass 
increase rate of the circumstellar disk is related to the initial stability of cloud core. We 
parameterized the ratio of thermal energy to gravitational energy, ao, which is related to the 
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mass accretion rate. In general, the mass accretion rate is described as 

M=fo.i, (24) 

where Mj is the Jeans mass and is the freefall timescale. In this study, we made an 
equihbrium sphere (i.e., critical BE sphere) and increased the density (or mass) by a factor 
of / (= 0.84/ao, see eq. [6]) to induce the collapse of the cloud core. Thus, the mass accretion 
rate in our setting can be described as 

where the initial mass of cloud core Md is proportional to the density enhancement factor 
/ and Jeans mass Mj as Md oc / Mj oc Mj, and the freefall timescale of the cloud core 
tff^ci is proportional to tfr_ci oc /~^'^^%,o oc a^^'^ts^, where ts^ is the freefall timescale of a 
critical BE sphere. Thus, a more unstable cloud core that has a smaller has a larger mass 
accretion rate and tends to show fragmentation due to gravitational instability. In reality, 
the mass increase rate for the model AOS is about three times higher than that for the model 
A05. Thus, the fragmentation and subsequent planets or substellar-mass objects may appear 
in relatively unstable cloud cores. 

In general, however, the initial clouds are characterized by three parameters, ao, Po, 
and 7o, as described in §2.2[ Although we fixed the angular velocity (/3o) and magnetic 
field strength (70) of the initial cloud core in this study, fragmentation may occur in more 
stable cloud cores (i.e ., cloud cores with large Oq) with different f3o and 70. As shown in 



Hennebelle fc Teyssierl (j2008bl ). in the collapsing cloud core, the magnetic field suppresses 



fragmentation while the rotation promot es it. Thus, fragment ation tends to occur in the 



cloud core with larger /3o and smaller 70 ( iMachida et al.ll2008al ). Moreover, even in a cloud 



core with a very weak magnetic field, fragmentation in the circumstellar disk may be sup- 
pressed by the global spiral structure that effectively transports the angular momentum and 
lowers the surface density (i.e., increases the Toomre Q-parameter) . On the other hand, 
when the circumstellar disk is strongly magnetized, the angular momentum is transferred 
by magnetic effects and the spiral structure hardly develops. Then, the mass accumulates 
in the inner disk region where Ohmic dissipation significantly weakens the magnetic field, 
and fragmentation can occur because no global spiral structure develops. Thus, magnetic 
field and rotation contribute to both promotion and suppression of fragmentation through a 
complicated process. In this study, we showed the possibility of fragmentation in the main 
accretion phase. However, to quantitatively determine fragmentation and subsequent forma- 
tion of planets and substellar-mass objects in the circumstellar disk, we should investigate 
the evolution of the circumstellar disk from the molecular cloud core in a large parameter 
space. 
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In the circumstellar disk, fragmentation occurs and planets or substellar-mass objects 
appear. In this study, we resolved the fragmentation process with sufficient spatial resolution 



( iTruelove et al. 1 119971 ). However, more spatial resolution and further developed numerical 
techniques may be necessary to determine the final mass of fragments. At its formation, 
a fragment's mass almost corresponds to the Jeans mass. In this study, we assumed the 
barotropic equation of state, in which the minimum Jeans mass is limited to ~ 5Mjup. 
Thus, fragmentation does not occur with a fragment mass of M ^ 5Mjup. If the circum- 
stellar disk cools to reaches lower temperatures, more less-massive fragments may appear. 
The barotropic approximation makes it possible to calculate very long-term evolution of 
the circumstellar disk and to investigate the evolution of the circumstellar disk by the end 
of the main accretion phase from the prestellar core stage, while radiation hydrodynamic 
calculations are necessary to qualitatively estimate the fragment mass. 

In addition, after fragmentation occurs in the circumstellar disk, we need a higher 
spatial resolution to investigate the acquisition process of the fragments' mass. Recent 
works about a gas giant planet formation showed that the Hill radius of the protoplanet 
should be resolved with sufficie nt spatial resolution for accurate estimate the rate of mass 



accretion onto the pr otoplanet (IKley et al. I l200ll : iD'Angelo et al. I l2003l : iBate et al.l 12003 



Machida et al.ll2010bl ). Note that in these calculations, the protoplanet mass was fixed, and 
the evolution of the circumstellar disk was calculated with a high spatial resolution for a 
shorter duration after a gas-giant formed. The Hill radius is expressed as 



where is the orbital radius of the planet. When the protoplanet orbiting the solar- mass 
protostar has a mass of Mp = 10~^Mps (10~^Mps) and is located at r = lAU (10 AU), the 
Hill radius is rn — 0.07 AU (1.5AU). When the Hill radius is resolved with ^ 10 cells, we 
need a minimum spatial resolution (i.e., the size of each cell) of ~ 0.007 AU (~ 0.15AU). 
However, calculation with such high resolution needs much computing power, and thus we 
only calculated the evolution of the circumstellar disk for a very short duration. In reality, 
the evolution of the circumstellar disk was calculated for only ~ 1000 yr at best in previous 
studies, while we calculated the evolution of the circumstellar disk from prestellar core stage 
until ~ 10^ yr after protostar formation. 

The final mass of fragments is important to determine whether they are planets or binary 
companions. To understand both formation and evolution of (gas-giant) planets, we need 
to calculate the evolution of the circumstellar disk for a longer duration with higher spatial 
resolution. However, such calculations are impossible even using present supercomputers 
because a longer temporal calculation confiicts with higher spatial resolution. Thus, in 
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this study, we covered the whole prestellar cloud core and investigated the evolution of 
the circumstellar disk with moderately coarser spatial resolution A ^ 0.6 AU, and showed 
that fragmentation occurs in the circumstellar disk and a planet or substellar-mass object 
forms. However, since we need a higher spatial resolution to determine the final mass of the 
fragment, we did not investigate this point in greater depth. In this study, we showed the 
possibility of the formation of planets or substellar-mass objects in a magnetized disk for the 
main accretion phase, and the precise mass of such objects should be investigated in future 
studies. 



4.3. Ohmic Dissipation and Ambipolar Diffusion 



As described in §3.6[ the recurrent planet formation and intermittent outflow shown 
in this paper are closely related to the effect of magnetic field and its dissipation in the 
circumstellar disk. Figure [TU] shows that the magnetic field is well-coupled with the neutral 
gas in the outer disk region {r^ > 10 — 100 AU) where a relatively high ionization rate is 
realized owing to a relatively low gas density. Thus, the angular momentum in such region 
is effectively transferred by the magnetic braking and protostellar outflow, and the gas can 
flow into the inner disk region. On the other hand, the magnetic field is not well-coupled 
with the neutral gas in the inner disk region (rc < 10 — 100 AU) where the gas density is high 
and a considerably low ionization rate is realized. Thus, the angular momentum transfer 
due to the magnetic effects is not so effective, and the gas flowing from the outer disk region 
accumulates in the inner disk region around the boundary between the magnetically active 
and inactive region at ~ 5 — 50 AU (see. Fig. [TOj) . In that region, the surface density 
continues to increase and fragmentation tends to occur. In summary, magnetic dissipation is 
essential in inducing recurrent fragmentation (or planet formation) in the circumstellar disk. 

In this study, we have considered only Ohmic dissipation as the dissipation mechanism 
of the magnetic field. This is because Ohmic dissipation is expected to be the most ef- 
fective mechanism for magnetic field dissipation in high density region where most of the 
magnetic flux leakage occurs. However, ambipolar diffusion may also contribute to the 
magnetic flux leakage in some intermediate density region. Recent studies have shown 
that the Ohmic dissipation is more effective than the ambipolar diffusion in the range of 
n ^ W 

Kunz fc Mouschoviad 120101 ) . In the range of 



10^2 cm"^ fiTassis fc MouschoviasI 12007 



n ^ 10 cm the magnetic field is well-coupled with the neutra l gas in a magnetically 



super critical cloud ( jNakano et al.ll2002l : iKunz &: MouschoviasI l2010l ). Thus, it might be im- 
portant to study how the inclusion of ambipolar diffusion would change the formation and 
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evolution of the circumstellar disk, through the possible contribution in the limited range of 
density. 

With strong coupling approximation, the induction equation including the ambipolar 
diffusion can be written as 



OB 

'dt 



V X (v X B) + tiodV^B + V X 



B 



X X (V X B)] 



(27) 



where p„, pi and 7 are the density of neutral gas and charged particles and a drag coe fficient 
between charged particles and neutrals, respectively (for details, see, e.g., IShulll992l ). The 
third term of equation f l27p can be regarded as a non-linear diffusion term by the ambipolar 
diffusion. When we adopt an approximate relation between neutral and charged particle 

1 /2 

density, pi = Cpn , where C is a constant, the effective diffusion coefficient of the ambipolar 
diffusion (rj^d) can be roughly estimated as 



'7ad 



47r7Cp3/2- 



(28) 



In a realistic setting, iKunz fc MouschoviasI ( I2OIOI ) showed that diffusion coefficient of the 
ambipolar diffusion (r^ad) is 10-100 times larger than that of the Ohmic dissipation (r^oo) 
in the range of 10^^ cm~^ ^ n <^ lO^^cm"^. Strictly speaking, to include the effect of the 
ambipolar diffusion, we should solve equation ( 1271) together with the equations determin- 
ing charged particle density. Note that the realistic evolution of charged particle density 
sensitively depends on total surface area of dust grains, since the recombination of charged 
particles on the grain surface is very efficient. However, there is a considerable uncertainty 
in the evolution of dust grain properties, such as the size distribution, in the protostellar 
collapse process. Thus, our limited knowledge of grain properties make our accurate treat- 
ment of magnetic field dissipation impracticable. This is also the reason why we adopted 
parametr ic description of the Ohmic diffusion coefficient in our study of magnetic field dis- 
sipation ( iMachida et al.l 120071 ). In this paper we take the same approach to the effect of 
ambipolar diffusion. To roughly estimate the effect of the ambipolar diffusion, we calculate 
the evolution of the circumstellar disk with an artificially larger resistivity. In equation ([3]), 
we replace r^oo by ?7od as follows: 



OB 

'of 



V X {v X B)+Cr, ?7odV^B. 



(29) 



We adopt = (ideal MHD model), 1 (fiducial model), 10, 100. 



As a first test, we constructed a non-rotating Bonner- Ebert sphere with a central density 
of n 10^ cm~^ immersed in the uniform magnetic field of Bz.o = 5 x 10~^ G. Then, we 
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calculated the collapse of this cloud with different c^. Figure [15] shows the evolution of the 
magnetic field as a function of the central number density. The magnetic field in each model 
tracks the same path in the range of n ^ 10^ cm~^, while the magnetic field for model with 
= 100 begins to depart from that for ideal MHD model (c^ = 0) at ra ~ lO^cm"^. The 
figure indicates that non-ideal MHD effect becomes important at the lower density with 
the larger c^. When the central number density reaches n = lO^^cm"^, the magnetic field 
strengths for model with = and 100 are = 7G and 0.01 G, respectively. Thus, the 
ideal MHD model has about 100 times stronger magnetic field than model with = 100 at 
this epoch. At the same epoch, models with = 1 and 10 have Bz = 0.03 G and 0.06 G of 
the magnetic field, respectively. Thus, the difference of the magnetic field strength between 
Cr] = 1 and 100 (10) is only factor of 6 (2) . 

Next, we calculated the evolution of the circumstellar disk with the initial condition 
same as model AOS but different resistivities = 0, 10 and 100. The density distribution 
on the equatorial plane after the circumstellar disk formation for models with = and 
Crj = 100 are plotted in Figure [161 As seen in this figure, fragmentation does not occur in 
ideal MHD model (c^ = 0) until the end of the calculation, while, for model with = 100, 
fragmentation occurs at — 2550 yr after the protostar formation and two clumps appear in 
the circumstellar disk. For model with Crj = 10, the circumstellar disk also shows fragmen- 
tation and subsequent planet formation. Note that we have terminated the calculation with 
larger resistivity earlier, since it takes formidable time to follow the evolution with large r 
resistivity with our time-explicit integration scheme for diffusion term (iMachida et al.ll2007l ). 
Thus, we cannot directly compare the results among models with different resistivities. How- 
ever, Figure [16] indicates that fragmentation occurs even when the magnetic field dissipates 
in a relatively low-density region. Therefore we expect that recurrent planet formation and 
intermittent outfiow can occur even with more realistic treatment of the ambipolar diffusion. 



4.4. Comparison with Previous Studies 



Only a few studies foc used on the fo r mation and evolution of the cir c umste llar disk from 
the prestellar core stage. IWalch et al.l ( l2009al ) and IVorobyov fc Basul (|2010[ ) investigated 
the circumstellar disk in unmagnetized clouds, and found that massive cir cumstellar disk 



formed at the early rn ain accretion phase tends to show fragmentation. IMachida et al. 
( I2OIOI ) and iBatd ( I2OIOI ) pointed out that the disk forms before the protostar formation, and 
the circumstellar disk is inevitably massive than the protostar in the early accretion phase. 
This is because the circumstellar disk is originated from the first ad i abatic core formed i n 
the gas-collapsing phase before the protostar formation ( iBatd Il998l ; IMachida et al.ll2010l ). 
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Such massive circumstellar disk is prone to show fragmentation. These studies indicate 
the possibihty of the formation of substellar- or planet-mass object in the unmagnetized 
circumstellar disk by the gravitational instability. 



In the magnetized cloud core, only I Vorobyov fc Basul (l2005l . l2006l ) investigated the for- 
mation and evolution of the circumstellar disk using two-dimensional simulation. Although, 
they adopted the magnetic evolution in an approximate way instead of solving the induction 
equation, they presented several interesting results. They pointed out that fragmentation 
and planet formation is possible even in the magnetized circumstellar disk. They also showed 
that the gravitational instability or fragments temporarily amplifies the mass accretion rate 
onto the protostar that may account for the FU Orionis outbursts. There are slight quantita- 
tive differences between I Vorobyov &: Basij (l2005l . l2006l ) and ours, because we calculated the 
circumstellar disk with a realistic process of the magnetic dissipation in three dimensions. 
However, our results are qualitatively the same as their results. This is because the circum- 
stellar disk alw a ys ex periences gravitationally unstable phase after the protostar formation. 
Inutsuka et al.l (120091 ) showed that, even in the magnetized cloud core, the circumstellar 
disk is massive than the protostar and becomes gravitationally unstable. Such disk tends to 
show fragmentation in the early accretion phase. Even when fragmentation does not occur 
in the early accretion phase, such massive disk is expected to show fragmentation in the 
later accretion or subsequent phase. Thus, we expect that the planet formation due to the 
gravitational instability is common in the star formation process. 



5. Summary 

To investigate the properties of the circumstellar disk during the main accretion phase, 
we calculated the formation and evolution of the circumstellar disk and protostellar outflow 
in a magnetized collapsing cloud core from the prestellar core stage using three-dimensional 
nested-grid simulation that covers both the whole region of a nascent molecular cloud core 
and sub-AU structure. We constructed two models with a parameter for the initial stability 
of cloud core (the ratio of thermal energy to gravitational energy) or total mass of the cloud 
core, in which the rotation rate and magnetic field strength were fixed, and calculated the 
evolution of the circumstellar disk about ~ lO'^yr after protostar formation. The different 
initial stabilities affect the rate of mass accretion onto the circumstellar disk; that is, a 
relatively massive disk appears with a high mass accretion rate in a relatively thermally 
unstable cloud core. We obtained the following results from these calculations: 

• Recurrent Planet Formation in the Circumstellar Disk 

The circumstellar disk formed in a relatively unstable cloud core is considerably massive 
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because of the high rate of mass accretion onto the disk and is in a highly gravitationally 
unstable state. Thus, fragmentation occurs and substellar-mass objects appear in the 
circumstellar disk. The substellar-mass objects lose their orbital angular momenta by 
gravitational interaction with the massive circumstellar disk and finally fall onto the 
protostar for ~ 1000 yr after their formation. After those objects fall, the circumstellar 
disk increases its mass by the mass accretion from the nascent cloud core and again 
shows fragmentation and subsequent falling of fragments onto the protostar. 

Fragmentation and falling are repeated in the main accretion phase. The mass of frag- 
ments never reaches the hydrogen-burning limit of ~ 0.08 Mq because the circumstellar 
disk self-regulates its mass during this recurrent planet formation process, that is, the 
circumstellar disk induces fragmentation and loses its mass when fragments fall onto 
the protostar before it becomes vary massive. Thus, fragments formed in the main 
accretion phase hardly exceed the planet or brown-dwarf mass. It is expected that this 
recurrent planet formation process continues as long as the disk accretes mass at the 
high rate. However, mass accretion onto the disk decreases with time. Thus, the final 
generation of fragments formed in the outer disk region just before the mass accretion 
stops (or significantly decreases) may survive to beco nie long-lived planets or brown- 



dwar f s with long orbital per iods, as recently observed ( jKalas et al.ll2008[ iMarois et al. 



20081 : iThalmann et al.l 120091 ) 



Highly Time- Variable Mass Accretion 

In the main accretion phase, both models (relatively stable and unstable models) show 
non-steady mass accretion onto the central protostar. In a relatively stable cloud core 
that shows no fragmentation, the rate of mass accretion onto the protostar fluctuates 
in the range of ~ 10^^ to ~ 10~^ Mq yr^^ with a period of ~ 100 yr. Thus, highly time- 
variable mass accretion is seen in this model. The magnitude of the mass accretion rate 
almost corresponds to the theoretical prediction. The cycle of variation of the mass 
accretion rate (or the accretion luminosity) corresponds to the orbital period of the 
location of clear spiral arms that effectively transfer the angular momentum outward 
and promote protostellar mass accretion. The spiral arms only appear in the inner 
disk region where a considerably weak field is realized owing to the effective Ohmic 
dissipation, and the angular momentum is not effectively transferred by the magnetic 
effect. As a result, the excess angular momentum causes the mass to accumulate in 
the inner disk region and increases the non-axisymmetric perturbation that promotes 
the mass accretion onto the protostar. 

Instead, in the outer disk region where the magnetic field is well-coupled with neutrals, 
the non-asymmetric perturbation hardly grows because the angular momentum is ef- 
fectively transferred by the magnetic effect. In summary, because the spiral structure 
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induces the time-variable protostellar accretion, the typical period of the variation of 
the accretion luminosity corresponds to the orbital extent of the magnetically inactive 
zone where the spiral patterns appear. Thus, the time variation or the cycle of the 
accretion luminosity is related to the magnetic dissipation process that is connected 
to the magnetic flux problem. Although we adopted a very realistic value of the mag- 
netic resistivity and ionization degree, we may also deduce them by observations of the 
luminosity variation of Class and I objects. In addition, the time variation of the 
accretion luminosity caused by spiral patterns may e xplain the luminos ity dispersion 
of Class and I objects over 2-3 orders of magnitude (lEvans et al.l 120091 ). 



Episodic Mass Accretion 

Episodic mass accretion is seen in a relatively unstable cloud core. Fragmentation 
occurs and substellar-mass objects repeatedly appear in the circumstellar disk formed 
in an unstable cloud core. The substellar-mass objects disturb the circumstellar disk 
and induce a strong episodic brightening of the protostar. The mass accretion rate is 
considerably low <^ lO^^MQyr"^ when the substellar-mass objects orbit in a region 
sufficiently far 5 AU) from the protostar, while the rate increases to ^ 10~^ Mq as 
the substellar-mass objects approach 5AU) the protostar. Thus, during the main 
accretion phase, the protostar shows two dis tinct phases of ma s s accr etion: quiescent 
and active mass accretion phases, as seen in IVorobyov &: Basd (120061 ). The protostar 
increases its luminosity {'^ 1 Lq) by a significant mass accretion in the active phase that 
lasts only for ~ 100 yr, while it dims 0.1 Lq) in the quiescent phase that lasts for 
;^ 1000 yr . Recent observations support episodic rn ass accretion to solve the luminosity 
problem (lEnoch et al.ll2009bl : iDunham et al.ll2010l ). In addition, FUori-type outbursts 
also may be related to the episodic accretion caused by orbital motion and falling of 
the substellar-mass objects. However, we need further spatial resolution to associate 
them with actual observations. 

Intermittent Protostellar Outflows 

Both models show the outfiow driven by the circumstellar disk. The steady outflow 
that lasts for > 10^ yr appears in a relatively stable cloud core, while the non-steady 
outflow appears in a relatively unstable cloud core. In a relatively stable cloud core, 
spiral arms that promote mass accretion onto the protostar appear in the inner disk 
region where magnetic dissipation is effective and no outflow appears, while the outflow 
is driven by the outer disk region where an almost axisymmetric structure remains. As 
a result, time- variable mass accretion is caused by the inner disk region, while outflow 
is steadily driven by the stable outer disk region. 

On the other hand, in a relatively unstable cloud core, the outflow driven by the cir- 
cumstellar disk shows a strong time variability because substellar-mass objects disturb 
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magnetic field lines in the whole circumstellar disk. As a result, the episodic (or in- 
termittent) outflow appears with episodic mass accretion. Thus, intermittent outflow 
may be proof of the existence of planets or substellar-mass objects in the circumstellar 
disk. In addition, variable periods of the outflow well synchronize with the orbital 
period of substellar-mass objects. Thus, we may be able to predict the existence of a 
planet and its orbital period with detailed observation of the trail of the protostellar 
outflow. However, further long-term calculations are necessary to investigate outflow 
propagation in the interstellar medium. 



Numerical computations were carried out on NEC SX-9 at Center for Computational 
Astrophysics, CfCA, of National Astronomical Observatory of Japan, and NEC SX-8 at 
the Yukawa Institute Computer Facility. This work was supported by Grants-in-Aid from 
MEXT (20540238, 21740136). 
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Table 1: Model parameters and calculation results 



Model 


QO 


/ 


ilo [s-l] 


Bini [mG] Mcl [Mq] 


Rc [AU] 


A 


Mps [Mq] 




Npi 


AOS 


0.3 


2.8 


1.1 X 10"^^ 


37 1.3 


4700 


9 


0.59 


0.15 


16 


AOS 


0.5 


1.7 


1.1 X 10-^3 


37 0.8 


4700 


5.6 


0.15 


0.03 






Table 2: Properties of planets 



Generation 


n [AU] 


ti [yr] 


te [yr] 


te — ti [yr] 


Mpi [Mjup] 


1 


6.8 


795 


3657 


2862 


23 


2 


45.8 


3990 


4669 


679 


6 


3 


22.5 


5037 


5693 


656 


26 


4 


36.9 


5984 


6603 


619 


16 


5 


36.5 


6697 


7503 


806 


15 


6 


40.2 


7504 


8446 


942 


9 


7 


23.3 


8455 


9240 


785 


27 


8 


25.9 


9327 


9760 


433 


25 


9 


22.1 


9748 






12 
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Fig. 1. — Time sequence images of the circumstellar disk after protostar formation for model 
A03. The density distribution on the equatorial plane around the protostar is plotted in each 
panel. The size of the black circle at the center in each panel corresponds to the sink radius. 
The box size of the upper three panels is ~ 80 AU, while that of the middle and lower panels 
is ~ 320 AU. The white squares in the middle and lower panels denote the outer boundary 
of the subgrid. The time elapsed after the protostar formation tc is plotted on the upper 
side of each panel. 



100 AU 



Fig. 2. — Configuration of protostellar outflow at tc = 843 yr is shown by yellow volume, in 
which color indicates outflow speed. The density distributions (colors) are projected on each 
wall surface. The velocity vectors (arrows) on the equatorial plane are plotted on bottom 
surface. The magnetic field lines integrated from each planet are plotted by yellow lines. 
The box size is 100 AU. 
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Fig. 3. — Density distribution of circumstellar disk on equatorial plane tc = 9755 yr after 
protostar formation. Eighth- and ninth-generation planets are shown. 
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Fig. 4. — Orbital evolution of the most massive planet against time elapsed after protostar 
formation. Symbol Q indicates the epoch at which fragmentation occurs and planet appears, 
while symbol x is the epoch at which planet falls onto protostar. Numbers 1-9 corresponds 
to generation of planet. 
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Fig. 5. — Trajectory for each planet just after planet formation (diamonds) until planet falls 
onto protostar (Ist-Sth-generation planets) or end of calculation (9th-generation planet). 
Size of black circle at center corresponds to sink radius. 
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Fig. 6. — Density distribution (color) on ?/ = plane (left panel) and z = plane (right 
panel) for model AOS at tc = 3781 yr after protostar formation. Arrows indicate the velocity. 
The left panel has a box size of 320 AU {I = 10, where / is the grid level), while the right 
panel has a box size of 160 AU (/ = 11). 
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Fig. 7. — Protostellar mass for models AOS (red) and AOS (black) against time elapsed after 
protostar formation. 
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Fig. 9. — Mass accretion rate for models AOS (black) and AOS (red) against time elapsed 
after protostar formation. 
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Fig. 10. — Radial distribution of the surface density (a), Toomre Q parameter (6), plasma 
beta (3p (c), and ratio of the azimuthal to the Kepler velocity (d) for models with AOS 
{tc = 445 and 3917yr) and AOS {tc = 3781 yr). 
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Fig. 11. — Outflowing mass for models A05 (black) and AOS (red) against time elapsed after 
protostar formation. Outflow is defined as gas with radial velocity exceeding sound speed 
Vr > Cs,o inside r < Re- 
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Fig. 12. — Ratio of outflowing to inflowing mass for models AOS (black) and AOS (red) 
against time elapsed after protostar formation. The outflowing mass is estimated as the 
mass ejected from the I = 9 grid, and the inflowing mass is estimated as mass entering the 
/ = 9 grid. 
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Fig. 13. — Accretion luminosity of protostar for models AOS (black) and AOS (red) against 
time elapsed time after protostar formation. 
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Fig. 14. — Mass accretion rate (left axis; black line and dots) and orbital radius of the most 
massive fragment (right axis; red line and dots) for model AOS during tc = 1000 — 2800 yr 
after protostar formation. The perihelion of the planet corresponds to the peak of the mass 
accretion rate. 
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Fig. 15. — The magnetic flux density against central number density for models with = 
(ideal MHD model), 1 , 10 and 100. 
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Fig. 16. — Density distribution and velocity vectors (arrows) on the equatorial plane for 
models = (ideal MHD; left panel) and 100 (right panel). The elapsed time after the 
protostar formation tc is described in the upper part of each panel. 



